Administrative boundaries regionizer

In [1]:

import geopandas as gpd
from osmnx import geocode_to_gdf
import plotly.express as px
from shapely.geometry import Point, box

from srai.regionizers import AdministrativeBoundaryRegionizer

import geopandas as gpd from osmnx import geocode_to_gdf import plotly.express as px from shapely.geometry import Point, box from srai.regionizers import AdministrativeBoundaryRegionizer

/opt/hostedtoolcache/Python/3.8.14/x64/lib/python3.8/site-packages/spherical_geometry/great_circle_arc.py:365: RuntimeWarning: invalid value encountered in divide
  return P / l

Regionize city¶

Basic usage of the AdministrativeBoundaryRegionizer with a city boundary.

Here admin_level equal to 9 defines city districts in Poland.

In [2]:

wroclaw_gdf = geocode_to_gdf(query=["R451516"], by_osmid=True)
wroclaw_gdf.plot()

wroclaw_gdf = geocode_to_gdf(query=["R451516"], by_osmid=True) wroclaw_gdf.plot()

Out[2]:

<AxesSubplot: >

In [3]:

abr = AdministrativeBoundaryRegionizer(admin_level=9)

abr = AdministrativeBoundaryRegionizer(admin_level=9)

In [4]:

wro_result_gdf = abr.transform(gdf=wroclaw_gdf)
wro_result_gdf.head()

wro_result_gdf = abr.transform(gdf=wroclaw_gdf) wro_result_gdf.head()

Loading boundaries: 0it [00:00, ?it/s][overpass] downloading data: [timeout:60][out:json];(relation["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="9"](51.0426686,16.8073393,51.2100604,17.1762192);); out ids tags;
Loading boundaries: 48it [00:55,  1.15s/it]

Out[4]:

	geometry
region_id
Osiedle Leśnica	POLYGON ((16.90228 51.18284, 16.90247 51.18202...
Osiedle Jerzmanowo-Jarnołtów-Strachowice-Osiniec	POLYGON ((16.83783 51.12350, 16.83855 51.12365...
Osiedle Pracze Odrzańskie	POLYGON ((16.91883 51.19454, 16.92015 51.19200...
Osiedle Żerniki	POLYGON ((16.93916 51.12235, 16.93851 51.12135...
Osiedle Maślice	POLYGON ((16.91968 51.15787, 16.92075 51.15993...

In [5]:

fig = px.choropleth_mapbox(
    wro_result_gdf,
    geojson=wro_result_gdf,
    color=wro_result_gdf.index,
    locations=wro_result_gdf.index,
    center={"lat": 51.125, "lon": 16.99},
    mapbox_style="carto-positron",
    zoom=10.5,
)
fig.update_layout(margin={"r": 0, "t": 0, "l": 0, "b": 0})
fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choroplethmapbox"))
fig.update_traces(showlegend=False)
minx, miny, maxx, maxy = wroclaw_gdf.geometry[0].bounds
fig.update_geos(
    projection_type="equirectangular",
    lataxis_range=[miny - 0.1, maxy + 0.1],
    lonaxis_range=[minx - 0.1, maxx + 0.1],
    showlakes=False,
    showcountries=False,
    showframe=False,
    resolution=50,
)
fig.update_layout(height=600, width=800, margin={"r": 0, "t": 0, "l": 0, "b": 0})
fig.show(renderer="png")  # replace with fig.show() to allow interactivity

fig = px.choropleth_mapbox( wro_result_gdf, geojson=wro_result_gdf, color=wro_result_gdf.index, locations=wro_result_gdf.index, center={"lat": 51.125, "lon": 16.99}, mapbox_style="carto-positron", zoom=10.5, ) fig.update_layout(margin={"r": 0, "t": 0, "l": 0, "b": 0}) fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choroplethmapbox")) fig.update_traces(showlegend=False) minx, miny, maxx, maxy = wroclaw_gdf.geometry[0].bounds fig.update_geos( projection_type="equirectangular", lataxis_range=[miny - 0.1, maxy + 0.1], lonaxis_range=[minx - 0.1, maxx + 0.1], showlakes=False, showcountries=False, showframe=False, resolution=50, ) fig.update_layout(height=600, width=800, margin={"r": 0, "t": 0, "l": 0, "b": 0}) fig.show(renderer="png") # replace with fig.show() to allow interactivity

Regionize country¶

How to return an empty region covering water bodies outside of the land.

Here admin_level equal to 4 defines country regions in Madagascar.

In [6]:

madagascar_bbox = box(minx=43.2541870461, miny=-25.6014344215, maxx=50.4765368996, maxy=-12.0405567359)
madagascar_bbox_gdf = gpd.GeoDataFrame({"geometry": [madagascar_bbox]}, crs="EPSG:4326")

madagascar_bbox = box(minx=43.2541870461, miny=-25.6014344215, maxx=50.4765368996, maxy=-12.0405567359) madagascar_bbox_gdf = gpd.GeoDataFrame({"geometry": [madagascar_bbox]}, crs="EPSG:4326")

In [7]:

abr = AdministrativeBoundaryRegionizer(admin_level=4, return_empty_region=True)

abr = AdministrativeBoundaryRegionizer(admin_level=4, return_empty_region=True)

In [8]:

madagascar_result_gdf = abr.transform(gdf=madagascar_bbox_gdf)
madagascar_result_gdf.tail()

madagascar_result_gdf = abr.transform(gdf=madagascar_bbox_gdf) madagascar_result_gdf.tail()

Loading boundaries: 0it [00:00, ?it/s][overpass] downloading data: [timeout:60][out:json];(relation["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="4"](-25.6014344215,43.2541870461,-12.0405567359,50.4765368996);); out ids tags;
Loading boundaries: 25it [00:35,  1.42s/it]

Out[8]:

	geometry
region_id
Atsinanana	POLYGON ((47.67041 -20.06128, 47.67343 -20.044...
Diana Region	MULTIPOLYGON (((48.35982 -13.60568, 48.35872 -...
Analanjirofo	MULTIPOLYGON (((48.90954 -16.67829, 48.92118 -...
Sava	POLYGON ((49.06377 -14.19881, 49.06448 -14.198...
EMPTY	MULTIPOLYGON (((50.47654 -25.60143, 45.21301 -...

In [9]:

fig = px.choropleth(
    madagascar_result_gdf,
    geojson=madagascar_result_gdf.geometry,
    locations=madagascar_result_gdf.index,
    color=madagascar_result_gdf.index,
    color_continuous_scale=px.colors.sequential.Viridis,
)

fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choropleth"))
fig.update_layout(coloraxis_showscale=False)
fig.update_traces(showlegend=False)
minx, miny, maxx, maxy = madagascar_bbox.bounds
fig.update_geos(
    projection_type="equirectangular",
    lataxis_range=[miny - 0.1, maxy + 0.1],
    lonaxis_range=[minx - 0.1, maxx + 0.1],
    visible=False,
    showlakes=False,
    showcountries=False,
    showframe=False,
    resolution=50,
)
fig.update_layout(height=800, width=420, margin={"r": 0, "t": 0, "l": 0, "b": 0})
fig.show(renderer="png")  # replace with fig.show() to allow interactivity

fig = px.choropleth( madagascar_result_gdf, geojson=madagascar_result_gdf.geometry, locations=madagascar_result_gdf.index, color=madagascar_result_gdf.index, color_continuous_scale=px.colors.sequential.Viridis, ) fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choropleth")) fig.update_layout(coloraxis_showscale=False) fig.update_traces(showlegend=False) minx, miny, maxx, maxy = madagascar_bbox.bounds fig.update_geos( projection_type="equirectangular", lataxis_range=[miny - 0.1, maxy + 0.1], lonaxis_range=[minx - 0.1, maxx + 0.1], visible=False, showlakes=False, showcountries=False, showframe=False, resolution=50, ) fig.update_layout(height=800, width=420, margin={"r": 0, "t": 0, "l": 0, "b": 0}) fig.show(renderer="png") # replace with fig.show() to allow interactivity

Regionize Europe¶

Option to slightly increase the value of toposiplify to simplify geometries even more.

Here admin_level equal to 2 defines countries.

In [10]:

eu_bbox = box(minx=-10.478556, miny=34.633284672291, maxx=34.597916, maxy=70.096054)
eu_bbox_gdf = gpd.GeoDataFrame({"geometry": [eu_bbox]}, crs="EPSG:4326")

eu_bbox = box(minx=-10.478556, miny=34.633284672291, maxx=34.597916, maxy=70.096054) eu_bbox_gdf = gpd.GeoDataFrame({"geometry": [eu_bbox]}, crs="EPSG:4326")

In [11]:

abr = AdministrativeBoundaryRegionizer(admin_level=2, toposimplify=0.0005)

abr = AdministrativeBoundaryRegionizer(admin_level=2, toposimplify=0.0005)

In [12]:

eu_result_gdf = abr.transform(gdf=eu_bbox_gdf)

eu_result_gdf = abr.transform(gdf=eu_bbox_gdf)

Loading boundaries: 0it [00:00, ?it/s][overpass] downloading data: [timeout:60][out:json];(relation["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="2"](34.633284672291,-10.478556,70.096054,34.597916);); out ids tags;
Loading boundaries: 55it [01:27,  1.58s/it]

In [13]:

eu_result_gdf.head()

eu_result_gdf.head()

Out[13]:

	geometry
region_id
Netherlands	MULTIPOLYGON (((3.08610 51.57045, 3.09529 51.5...
Portugal	POLYGON ((-9.05070 41.86519, -8.87097 41.86393...
Morocco	POLYGON ((-5.88545 35.96488, -5.66652 35.93567...
Ireland	POLYGON ((-6.71943 55.28017, -6.82172 55.21647...
Faroe Islands	POLYGON ((-8.11704 62.11484, -8.11273 62.13026...

In [14]:

fig = px.choropleth(
    eu_result_gdf,
    geojson=eu_result_gdf.geometry,
    locations=eu_result_gdf.index,
    color=eu_result_gdf.index,
    color_continuous_scale=px.colors.sequential.Viridis,
)

fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choropleth"))
fig.update_layout(coloraxis_showscale=False)
fig.update_traces(showlegend=False)
minx, miny, maxx, maxy = eu_bbox.bounds
fig.update_geos(
    projection_type="equirectangular",
    lataxis_range=[miny - 1, maxy + 1],
    lonaxis_range=[minx - 1, maxx + 1],
    showlakes=False, 
    showcountries=False,
    showframe=False,
    resolution=50,
)
fig.update_layout(height=800, width=1000, margin={"r": 0, "t": 0, "l": 0, "b": 0})
fig.show(renderer="png")  # replace with fig.show() to allow interactivity

fig = px.choropleth( eu_result_gdf, geojson=eu_result_gdf.geometry, locations=eu_result_gdf.index, color=eu_result_gdf.index, color_continuous_scale=px.colors.sequential.Viridis, ) fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choropleth")) fig.update_layout(coloraxis_showscale=False) fig.update_traces(showlegend=False) minx, miny, maxx, maxy = eu_bbox.bounds fig.update_geos( projection_type="equirectangular", lataxis_range=[miny - 1, maxy + 1], lonaxis_range=[minx - 1, maxx + 1], showlakes=False, showcountries=False, showframe=False, resolution=50, ) fig.update_layout(height=800, width=1000, margin={"r": 0, "t": 0, "l": 0, "b": 0}) fig.show(renderer="png") # replace with fig.show() to allow interactivity

Toposimplify differences¶

Shows differences in simplification of small regions using four values: 1e-4, 1e-3, 1e-2 and 0.1. Those values are in degress, since it uses Douglas-Peucker simplification algorithm.

1e-4 is the default value and is equal to about 11.1m accuracy.

More info: https://github.com/mattijn/topojson

Here admin_level equal to 6 defines city districts in Singapore.

In [15]:

singapore_bbox = box(minx=103.5111238, miny=1.1263707, maxx=104.1313374, maxy=1.4787511)
singapore_bbox_gdf = gpd.GeoDataFrame({"geometry": [singapore_bbox]}, crs="EPSG:4326")

singapore_bbox = box(minx=103.5111238, miny=1.1263707, maxx=104.1313374, maxy=1.4787511) singapore_bbox_gdf = gpd.GeoDataFrame({"geometry": [singapore_bbox]}, crs="EPSG:4326")

In [16]:

results = {}
for value in [0.0001, 0.001, 0.01, 0.1]:
    abr = AdministrativeBoundaryRegionizer(admin_level=6, toposimplify=value)
    results[value] = abr.transform(gdf=singapore_bbox_gdf)

results = {} for value in [0.0001, 0.001, 0.01, 0.1]: abr = AdministrativeBoundaryRegionizer(admin_level=6, toposimplify=value) results[value] = abr.transform(gdf=singapore_bbox_gdf)

Loading boundaries: 0it [00:00, ?it/s][overpass] downloading data: [timeout:60][out:json];(relation["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="6"](1.1263707,103.5111238,1.4787511,104.1313374);); out ids tags;
Loading boundaries: 8it [00:09,  1.23s/it]
Loading boundaries: 8it [00:00, 88.34it/s]
Loading boundaries: 8it [00:00, 88.84it/s]
Loading boundaries: 8it [00:00, 77.28it/s]

In [17]:

minx, miny, maxx, maxy = singapore_bbox.bounds
for epsilon, result in results.items():
    fig = px.choropleth_mapbox(
        result,
        geojson=result,
        color=result.index,
        locations=result.index,
        center={"lat": 1.3119350704252704, "lon": 103.82412242562575},
        mapbox_style="carto-positron",
        zoom=9.5,
    )
    fig.update_layout(margin={"r": 0, "t": 0, "l": 0, "b": 0})
    fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choroplethmapbox"))
    fig.update_traces(showlegend=False)
    fig.update_geos(
        projection_type="equirectangular",
        lataxis_range=[miny - 0.1, maxy + 0.1],
        lonaxis_range=[minx - 0.1, maxx + 0.1],
        showlakes=False,
        showcountries=False,
        showframe=False,
        resolution=50,
    )
    
    size = len(result.to_json().encode('utf-8'))
    fig.update_layout(
        height=450,
        width=700,
        margin={"r": 0, "t": 50, "l": 0, "b": 0},
        title_text=f"Toposimplify value: {epsilon} ({size/1000} KB)",
    )

    fig.show(renderer="png")  # replace with fig.show() to allow interactivity

minx, miny, maxx, maxy = singapore_bbox.bounds for epsilon, result in results.items(): fig = px.choropleth_mapbox( result, geojson=result, color=result.index, locations=result.index, center={"lat": 1.3119350704252704, "lon": 103.82412242562575}, mapbox_style="carto-positron", zoom=9.5, ) fig.update_layout(margin={"r": 0, "t": 0, "l": 0, "b": 0}) fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choroplethmapbox")) fig.update_traces(showlegend=False) fig.update_geos( projection_type="equirectangular", lataxis_range=[miny - 0.1, maxy + 0.1], lonaxis_range=[minx - 0.1, maxx + 0.1], showlakes=False, showcountries=False, showframe=False, resolution=50, ) size = len(result.to_json().encode('utf-8')) fig.update_layout( height=450, width=700, margin={"r": 0, "t": 50, "l": 0, "b": 0}, title_text=f"Toposimplify value: {epsilon} ({size/1000} KB)", ) fig.show(renderer="png") # replace with fig.show() to allow interactivity

Regionize points¶

How to return original regions without clipping and select them using list of points. Showed using list of metro stations in Paris.

Here admin_level equal to 8 defines communes in France.

In [18]:

import requests

r = requests.get("https://raw.githubusercontent.com/w8r/paris-metro-graph/master/metro.json").json()
stations_gdf = gpd.GeoDataFrame({"geometry": [Point(s['longitude'], s['latitude']) for s in r['nodes']]}, crs="EPSG:4326")
stations_gdf

import requests r = requests.get("https://raw.githubusercontent.com/w8r/paris-metro-graph/master/metro.json").json() stations_gdf = gpd.GeoDataFrame({"geometry": [Point(s['longitude'], s['latitude']) for s in r['nodes']]}, crs="EPSG:4326") stations_gdf

Out[18]:

	geometry
0	POINT (2.30159 48.81081)
1	POINT (2.29703 48.81534)
2	POINT (2.24947 48.88824)
3	POINT (2.27210 48.88106)
4	POINT (2.28937 48.87558)
...	...
309	POINT (2.33812 48.88228)
310	POINT (2.34955 48.88729)
311	POINT (2.35482 48.84599)
312	POINT (2.29583 48.87493)
313	POINT (2.40655 48.87706)

314 rows × 1 columns

In [19]:

abr = AdministrativeBoundaryRegionizer(admin_level=8, return_empty_region=False, clip_regions=False)
paris_districts_result = abr.transform(gdf=stations_gdf)

abr = AdministrativeBoundaryRegionizer(admin_level=8, return_empty_region=False, clip_regions=False) paris_districts_result = abr.transform(gdf=stations_gdf)

Loading boundaries: 0it [00:00, ?it/s][overpass] downloading data: [timeout:60][out:json];is_in(48.8108097,2.3015874);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 1it [00:02,  2.07s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8153367,2.2970255);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 2it [00:04,  2.01s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8882396,2.2494747);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 3it [00:06,  2.05s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8810615,2.2720961);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 4it [00:08,  2.03s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8755816,2.2893655);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 5it [00:10,  2.05s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8336468,2.2435931);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 6it [00:12,  2.04s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8241474,2.2731713);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 7it [00:14,  2.06s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8445799,2.438912);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 8it [00:16,  2.05s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.846435,2.4184422);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 9it [00:18,  2.06s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8653357,2.416682);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 10it [00:20,  2.09s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.9163969,2.2946541);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 11it [00:22,  2.13s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8920499,2.2850318);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 12it [00:25,  2.16s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.857946,2.4357587);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 13it [00:27,  2.17s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.9032278,2.3059879);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 14it [00:29,  2.21s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8264064,2.406201);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 15it [00:31,  2.21s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.89128,2.4025383);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 16it [00:34,  2.31s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8797602,2.4164717);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 17it [00:36,  2.28s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.906204,2.331871);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 18it [00:38,  2.23s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.9231172,2.2862019);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 19it [00:40,  2.19s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.9142573,2.4038072);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 20it [00:42,  2.14s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.9207587,2.4107677);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 21it [00:44,  2.10s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.7963514,2.449235);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 22it [00:46,  2.08s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8033311,2.4456085);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 23it [00:48,  2.11s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.9195753,2.3433012);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 24it [00:50,  2.08s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8051317,2.3637903);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 25it [00:52,  2.06s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.810328,2.3622359);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 26it [00:54,  2.04s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.8159011,2.3773455);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 27it [00:57,  2.13s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.9064052,2.4491919);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 28it [00:59,  2.15s/it][overpass] downloading data: [timeout:60][out:json];is_in(48.818234,2.3195203);area._["boundary"="administrative"]["type"~"boundary|multipolygon"]["admin_level"="8"];relation(pivot); out ids tags;
Loading boundaries: 29it [01:02,  2.15s/it]

In [20]:

paris_districts_result.head()

paris_districts_result.head()

Out[20]:

	geometry
region_id
Châtillon	POLYGON ((2.27988 48.81091, 2.28325 48.81056, ...
Malakoff	POLYGON ((2.30132 48.82513, 2.31413 48.82226, ...
Puteaux	POLYGON ((2.25372 48.88691, 2.25244 48.88551, ...
Neuilly-sur-Seine	POLYGON ((2.24562 48.87636, 2.24777 48.87849, ...
Paris	POLYGON ((2.23208 48.86951, 2.24046 48.87189, ...

In [21]:

# Paris
fig = px.choropleth_mapbox(
    paris_districts_result,
    geojson=paris_districts_result,
    color=paris_districts_result.index,
    locations=paris_districts_result.index,
    center={"lat": 48.858, "lon": 2.353},
    mapbox_style="carto-positron",
    zoom=10.8,
)
fig2 = px.scatter_mapbox(
    stations_gdf, lat=stations_gdf.geometry.y, lon=stations_gdf.geometry.x
)
fig.add_trace(fig2.data[0])
fig.update_layout(margin={"r": 0, "t": 0, "l": 0, "b": 0})
fig.update_traces(showlegend=False)
fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choroplethmapbox"))
fig.update_traces(marker_color="black", marker_size=5, selector=dict(type="scattermapbox"))
fig.update_layout(height=800, width=800, margin={"r": 0, "t": 0, "l": 0, "b": 0})
fig.show(renderer="png")  # replace with fig.show() to allow interactivity

# Paris fig = px.choropleth_mapbox( paris_districts_result, geojson=paris_districts_result, color=paris_districts_result.index, locations=paris_districts_result.index, center={"lat": 48.858, "lon": 2.353}, mapbox_style="carto-positron", zoom=10.8, ) fig2 = px.scatter_mapbox( stations_gdf, lat=stations_gdf.geometry.y, lon=stations_gdf.geometry.x ) fig.add_trace(fig2.data[0]) fig.update_layout(margin={"r": 0, "t": 0, "l": 0, "b": 0}) fig.update_traces(showlegend=False) fig.update_traces(marker={"opacity": 0.6}, selector=dict(type="choroplethmapbox")) fig.update_traces(marker_color="black", marker_size=5, selector=dict(type="scattermapbox")) fig.update_layout(height=800, width=800, margin={"r": 0, "t": 0, "l": 0, "b": 0}) fig.show(renderer="png") # replace with fig.show() to allow interactivity